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ABSTRACT 

We compute the effect of primordial non-Gaussianity on the halo mass function, using excursion set theory. 
In the presence of non-Gaussianity the stochastic evolution of the smoothed density field, as a function of the 
smoothing scale, is non-markovian and beside "local" terms that generalize Press-Schechter (PS) theory, there 
are also "memory" terms, whose effect on the mass function can be computed using the formalism developed 
in the first paper of this series. We find that, when computing the effect of the three-point correlator on the 
mass function, a PS-like approach which consists in neglecting the cloud-in-cloud problem and in multiplying 
the final result by a fudge factor ~ 2, is in principle not justified. When computed correctly in the framework 
of excursion set theory, in fact, the "local" contribution vanishes (for all odd-point correlators the contribution 
of the image gaussian cancels the Press-Schechter contribution rather than adding up), and the result comes 
entirely from non-trivial memory terms which are absent in PS theory. However it turns out that, in the limit of 
large halo masses, where the effect of non-Gaussianity is more relevant, these memory terms give a contribution 
which is the the same as that computed naively with PS theory, plus subleading terms depending on derivatives 
of the three-point correlator. We finally combine these results with the diffusive barrier model developed in the 
second paper of this series, and we find that the resulting mass function reproduces recent A^-body simulations 
with non-Gaussian initial conditions, without the introduction of any ad hoc parameter. 

Subject headings: cosmology:theory — dark mattenhalos — large scale structure of the universe 



1. INTRODUCTION 
In th e firs t two papers of thi s series (Maggior e & Riottol 



d2009al) and iMaggiore & Riottol d2009bl) . papers I and II in 
the following) we have studied the mass function of dark mat- 
ter halos using the excursion set formalism. The halo mass 
function can be written as 



dn(M) „ s p d In a~HM) 



dM 



M 2 dlnM 



(D 



where n(M) is the number density of dark matter halos of 
mass M, a(M) is the variance of the linear density field 
smoothed on a scale R corresponding to a mass M, and p 
is the average density of the universe. Analytical deriva- 
tions of the halo mass fu nction are typically bas ed on Press- 
Schecht er (PS) theory dPress & Schechterlll974l) and its ex- 
tension dPeacock & Heavenslll990t iBond et al.ll 1 99 lb known 
as excursion set theory (see IZentnerl (120071) for a recent re- 
view). In excursion set theory the density perturbation evolves 
stochastically with the smoothing scale, and the problem of 
computing the probability of halo formation is mapped into 
the so-called first-passage time problem in the presence of a 
barrier. With this method, for gaussian fluctuations one ob- 



tains 



2 \ 1/2 <5, 



/PS(<7)= - - 

n J a 



■c -S 2 /(2a 2 ) 



(2) 



where 5 C — 1 .686 is the critical value in the spherical collapse 
model. This result can be extended to arbitrary redshift z by 
reabsorbing the evolution of the variance into S c , so that 8 C in 
the above result is replaced by S c (z) = 5 C (Q)/D(z), where D(z) 
is the linear growth factor. Equation (f2]i is only valid when the 
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density contrast is smoothed with a sharp filter in momentum 
space. In this case the evolution of the density contrast S(R) 
with the smoothing scale is markovian, and the probability 
that the density contrast reaches a given value S at a given 
smoothing scale satisfies a Fokker-Planck equation with an 
"absorbing barrier" boundary condition. From the solution 
of this equation one obtains eq. (0, including a well-known 
factor of two that Press and Schecther were forced to add by 
hand. 

However, as is well-known, a sharp filter in momentum 
space is not appropriate for comparison with experimental 
data from upcoming galaxy surveys, nor with A^-body simu- 
lations, because it is not possible to associate unambiguously 
a mass M to the smoothing scale R used in this filter. Rather, 
one should use a tophat filter in coordinate space, in which 
case the mass associated to a smoothing scale R is trivially 
(4/3)irR i p. If one wants to compute the halo mass func- 
tion with a tophat filter in coordinates space one is confronted 
with a much more difficult problem, where the ev olution of 6 
with t he smoothing scale is no longer markovian (Bo nd et ail 
1991). Nevertheless, in paper I we succeeded in developing a 
formalism that allows us to compute perturbatively these non- 
markovian effects and, for gaussian fluctuations, we found 
that, to first order, eq. (f2j) is modified to 



where 
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k(R)= lim 



(S(R')S(R)) 



1 -0.4592 -0.0031/?, 



(3) 



(4) 



W)> 

R is measured in Mpc/h, T(0,z) is the incomplete Gamma 
function, and the numerical value of k(R) is computed using a 
tophat filter function in coordinate space and a ACDM model 
with cr 8 = 0.8, h = 0.7, fl M = 1 - fi A = 0.28, n B = 0.046 and 
n s = 0.96, consistent with the WMAP 5-years data release. 
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This analytical result reproduces well the result of a Monte 
Carlo realization of the first-crossing distribution of excursion 
set theory, obtained by integrating numerically a Langevin 
equati on with a colored noise, per formed in iBond et all 
(119911) and in iRobertson et al.l (120081) . This is a useful test 
of our technique. Still, neither eq. (ffj) nor eq. ^ perform 
well when compared to cosmological A^-body simulation, 
which means that some crucial physical ingredient is still 
missing in the model. This is not surprising, since the 
spherical (or ellipsoidal) collapse model is a very simplified 
description of the process of halo formation which, as 
shown by N-body simulations, is much more complicated, 
and proceeds through a mixture of smooth accretion and 
violent encoun ters leading to merg ing as well as to frag- 
mentation (see Springel et al. (2005) and the related movies 
at http://www.mpa-garching.mpg.de/galform/millennium/). 
Furthermore, the very operative definition of what is a dark 
matter halo is a subtle issue. Real halos are not spherical. 
They are at best triaxial, and often much more irregular, expe- 
cially if they experienced recent mergers. Searching for them 
with a spherical overdensity algorithm therefore introduces 
further statistical uncertainties. Similar considerations hold 
for Friends-of-Friends halo finders. 

In paper II we have discussed in detail these uncertainties 
and we have proposed that at least some of the complications 
intrinsic to a realistic process of halo formation (as well as the 
statistical uncertainties related to the details of the halo finder) 
can be accounted for, within the excursion set framework, by 
treating the critical value for collapse as a stochastic variable. 
In this approach all our ignorance on the details of halo for- 
mation is buried into the variance of the fluctuations of the 
collapse barrier. The computation of the halo mass function 
is then mapped into a first-passage time process in the pres- 
ence of a diffusing barrier, i.e. a barrier whose height evolves 
according to a diffusion equation. For gaussian fluctuations 
we found that eq. <[3j must be replaced by 
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and Db is an effective diffusion coefficient for the barrier. A 
first-principle computation of Db appears difficult, but from 
recent studies of t he properties of the c ollapse barrier in N- 
body simulations (Robertson et al. 2008) we deduced a value 
D B — (0.3S C ) 2 . Using this value for D B in eq. © gives a ~ 
0.80, so 

^-0.89. (7) 

We see that the net effect of the diffusing barrier is that, in 
the mass function, S c is replaced by a x l 2 8 c , which is the re- 
placement that was ma d e by hand, simply to fit the data, in 
ISheth & Tormenl (fT999h : ISheth et al.l d200lh . 

The above result was obtained by considering a barrier that 
fluctuates over the constant value 5 C of the spherical collapse 
model. More generally, one should consider fluctuations over 
the barrier B{a) given by the ellipsoidal collapse model. Since 
the latter reduces to the former in the small a limit (i.e. for 
large halo masses), eq. (0 is better seen as the large mass 



limit of a more accurate mass function obtained from a bar- 
rier that flucutates around the average value B(a) given by 
the ellipsoidal collapse model. When k = eq. Q reduces to 
the large mass limit of the Sheth-Tormen mass function. So, 
eq. (O generalizes the Sheth-Tormen mass function by taking 
into account the effect of the tophat filter in coordinate space, 
while eq. © provides a physical motivation for the introduc- 
tion of the parameter a. 

Equation (0 is in excellent agreement with the iV-body sim- 
ulations for gaussian primordial fluctuations, see Figs. 6 and 
7 of paper II. We stress that our value a ~ 0.80 is not de- 
termined by fitting the mass function to the data. We do 
have an input from the Af-body simulation here, which is 
however quite indirect, and is the measured variance of the 
threshold for collapse , which for small a is determined in 
IRobertson et al.l (120081) to be ~ 0.3cr. Our diffusing barrier 
model of paper II translates this information into an effective 
diffusion coefficient for the barrier, Db = (0.3S C ) 2 , and predicts 
a = 1/(1 +Db). We refer the reader to paper II for details and 
discussions of the physical motivations for the introduction of 
a stochastic barrier. 

The above results refer to initial density fluctuations which 
have a gaussian distribution. In this paper we attack the 
problem of the effect on the halo mass function of non- 
Gaussianities in the primordial density field. Over the 
last decade a great deal of evidence has been accumulated 
from the Cosmic Microwave Background (CMB) anisotropy 
and Large Scale Structure (LSS) spectra that the observed 
structures originated from seed fluctuations generated dur- 
ing a primordial stage of inflation. While standard one- 
single field models of slow-roll inflation p redict that these 
fluctua t ions are very clos e to gaussian (see lAcquaviva et alj 
(HH; Maldacena (2003)), non-standard scenarios allow for 
a larger level of non-Gaussianity (see B artolo etail (120041) and 
refs. therein). Deviations from non-Gaussianity are usually 
param etrized by a dimensionless quantity /nl dBartolo et alJ 
(2004)) whose value sets the magnitude of the three-point 
correlation function. If the process generating the primor- 
dial non-Gaussianity is local in space, the parameter /nl in 
Fourier space is independent from the momenta entering the 
three-point correlation function; if instead the process is non- 
local in space, like in models of inflation with non-canonical 
kinetic terms, /nl acquires a dependence on the momenta. It 
is clear that detecting a significant amount of non-Gaussianity 
and its shape either from the CMB or from the LSS offers the 
possibility of opening a window into the dynamics of the uni- 
verse during the very first stages of its evolution. Current lim- 
its on the strength of no n-Gaussianity set the f m parameter to 
be smaller than 0(100) dKomatsu et alJd2008l) ). 

Non-Gaussianities are particularly relevant in the high- 
mass end of the power spectrum of perturbations, i.e. 
on the scale of galaxy clusters, since the effect of non- 
Gaussian fluctuations becomes especially visible on the 
tail of the probability distribution. As a result, both 
the abundance and the clustering properties of very 
massive halo s are sensitive probes of primordial non- 
Gaussianities (Matarrese et al. 1986; Grinstein & Wise 1986; 
Lucc hinet alJ 119881: iMo scardini e t al.l 11991b iKovama et . 
19991: iMatarrese et all 120001; iRobinson & Baker I |20Qi 



Robin son et 

aDHOOO), 

and could be detected or significantly 



constrained by the various planned large-scale galaxy sur- 
veys, both ground based (such as DES, PanSTARRS and 
LSS T) and on s a tellite (suc h as EUCLID and AD EPT) see, 
e.g. iDalal et alJ (120081) and ICarbone et ail (120081) . Further- 
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more, the primordial non-Gaussianity alters the clustering 
of dark matte r halos inducing a scale-dependent bias on 
large scales dDalal et all 120081; IMatarrese & Verdd 120081; 
ISlosar et alj l2008t lAfshordi & Tollevl 120081) whileeven for 
small primordial non-Gaussianity the evolution of perturba- 
tions on super -Hubble scales yiel ds extra contributions on 
smaller scales (Ba rtolo et al.l ( 120051) ). 

At present, there exist already various iV-body simulations 
where non-Gaussianity has been included in the initial condi- 
tions dKang et alj 120071: I Grossi et alj|2007t iDalal et alj|2008l: 
Desia cques et alj 120081: iPillepich et alj 120081; IGrossi et alj 
2009) and which are useful to test the accuracy of the different 
theoretical predictions for the dark matter halo mass function 
with non-Gaussianity. 

Various attempts at computing analytically the effect 
of primordial non-Gaussianities on the mass function ex- 
ist in the liter ature, based on non-Gaussian extensions 
of PS theory dChiu et all fl997l: iRobinson & Baker 1 120001 
iMatarrese etail l2000t lLo Verde et alJ 120081) . However, for 
gaussian fiucutations, in the large mass regime PS theory is off 
by one order of magnitude. It is clear that, by computing non- 
Gaussian corrections over a theory that, already at the gaus- 
sian level, in the relevant regime is off by an order of magni- 
tude, one cannot hope to get the correct mass function for the 
non-Gaussian case. What is typically done in the recent litera- 
ture is to take the ratio Rng(M) of the non-Gaussian halo mass 
function to the gaussian halo mass function, both computed 
within the framework of PS theory, hoping that even if nei- 
ther the former nor the latter are correct, still their ratio might 
catch the main modifications due to non-Gaussianities. The 
full non-Gaussian halo mass function is then obtained by tak- 
ing a fit to the data in the gaussian case, such as the Sheth and 
Thormen mass function (Sheth & Tormen 1999; Sheth et al. 
120011) . and multiplying it by R N g(M ). With this philosophy, 
the result of IMatarrese et all d2000l) reads 4 



^ng(ct) = exp 



S 3 c S 3 (a) 
6a 2 



(8) 



where 



6 ^/l-5 c S 3 (cr)/3dln<T 

<?M {S3(S)) 



+ y/\-S c S 3 (cT)/3 



(9) 



(S 2 (S)) 2 

is the (normalized) skewness of the density field and, as usual, 
S = a 2 is the variance. Since a = <r(M), we can equivalently 
consider /?ng as a function of M. 5 

With a similar philosophy, but a differ ent expansion 
techni que, namely the Edgeworth expansion, lLo Verde et al.1 
(2008) propose 

2 r 

T 4 
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In the limit a/S c <C 1, eq. (TTOb becomes 

S 3 c S 3 (a) 



Rng(<j) = 1 + ■ 



6a 2 



(10) 



(11) 



4 W e thank S. Matarrese for pointing out to us a typo in Matarr ese et alj 
<200Ch . 

5 We do not write explicitly the dependence of R^q(ij) on redshift z, which 
enters through the variance a 2 and, as usual, can be reabsorbed into the height 
6c of the critical val ue for collapse. Th e normalized skewness must instead 
be taken at z = 0, see IGrossi et all f2009). 



The same result is obtained from eq. <[8J expanding first to 
linear order in S 3 ((j), and then retaining the leading term of 
the expansion for small <r/8 c . The two formulas differ instead 
at the lev el of the term s suble ading in the expansion for small 
a/5 c . In IGrossi et ail d2009l) . in order to fit the data of N- 
body simulations, it was suggested to modify both eq. ([S]) and 
eq. ([Tol l, by making the replacement 

6 c ->5 e £ = ^/a6 c , (12) 
with a value ^[a ~ 0.86 obtained from the fit to the data, very 
close to our prediction given in eq. (|7]). 6 In the gaussian case 
we have shown in paper II that this replacement, which in 
the previous literature was made ad hoc just to fit the data, 
actually follows from the diffusive barrier model, see eq. ©, 
and that the precise value of a depends, among other things, 
on the details of the halo finder in the simulation, so (slightly) 
different values of a are obtained from A^-body simulations 
with different halo finders. Below we will see how the results 
of pa per II generalize to the non-Gaussian case. 

In lGrossi et al.1 d2009l) it is shown that, after performing the 
replacement dT2l) . both eq. (O and eq. dTOb are in good agree- 
ment with the result of A^-body simulations with non-Gaussian 
initial conditions, which a posteriori can be seen as a justifi- 
cation of the procedure used in their derivation. However, it is 
clear that taking the ratio of two results that, in the interesting 
mass range, are known to be both off by one order of mag- 
nitude, in order to get a fine effect such as the non-Gaussian 
corrections, can only be considered as a heuristic procedure. 
First of all, PS theory by itself produces a wrong exponential 
factor, since it would give a = 1 . Here one might argue that 
the gaussian and non-Gaussian mass functions have the same 
exponential behavior, so this effect cancels when considering 
the ratio R^g, and is anyhow accounted for by the heuristic 
prescription dT2l >. Still, a further source of concern is tha t 
the derivation of the PS mass function in iBond et ail d!991l) 
requires that the density field S evolves with the smoothing 
scale R (or more precisely with S(R)) in a markovian way. 
Only under this assumption one can derive eq. (0 together 
with the correct factor of two that Press and Schechter were 
forced to introduce by hand. As we have discussed at length 
in paper I, this markovian assumption is broken by the use 
of a filter function different from a sharp filter in momentum 
space and, of course, it is further violated by the inclusion of 
non-Gaussian corrections. When studying non-Gaussianities, 
it is therefore particularly important to perform the compu- 
tation including the effect of the filter, otherwise one would 
attribute to primordial non-Gaussianities effects on the mass 
functions which are due, more trivially, to the filter. 

The formalism that we have developed in papers I and II, 
however, allows us to attack the problem. First of all, in pa- 
per I we have set up a "microscopic" approach which is in 
principle exact. 7 With this formalism, we computed the non- 
markovian corrections due to the filter function, which are 

6 As discussed in paper II, the value of the diffusion constant of the barrier 
Db, and hence our prediction for a, depends on the halo finder. The value 
D b — 0.25, which leads to ^fa ~ 0.89, has been deduced from the simulation 
of lRobert son et al. 12008), that uses a s pherical overdensity (SO) halo finder 
with A = 200, while Grossi et al. 1 2009) use a friends-of-friends (FOF) halo 
finder with link-length 0.2. In the gaussian case, the mass functions obtained 
from these two finders are very close to each other. However, in order to 
perform an accurate numerical comparison of our prediction with TV-body 
simulations with non-Gaussian initial conditions, in would be necessary to 
determine both Dg and the mass function with the same halo finder. 

7 By exact we mean that, given the problem of halo formation as it is for- 
mulated within excursion set theory, the path integral technique developed in 
paper I is an exact way of attacking the mathematical problem of the first- 
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given by the terms proportional to k in eq. (0). This is impor- 
tant because it allows us to subtract, from a measurement of 
the halo mass function, the "trivial" effect due to the filter, and 
to remain with the effects due to genuine non-Gaussianities. 
Second, putting together the corrections due to the filter with 
the model of a diffusing barrier, we ended up with a halo 
mass function which works very well in the gaussian case, see 
Figs. 6 and 7 of paper II, and which therefore is a meaningful 
starting point for the inclusion of non-Gaussian perturbations. 
Finally, the formalism developed in paper I can be applied, 
with simple modifications, to the perturbative computation of 
the non-Gaussian corrections. This will be the subject of the 
present paper. 

This paper is organized as follows. In Section|2j extending 
to the non-Gaussian case the results presented in paper I, we 
show how to formulate the first-passage time problem for non- 
Gaussian fluctuations in terms of a path integral with bound- 
aries, and we recall the basic points of the computation of non- 
markovian corrections performed in paper I. In Section [3] we 
compute the non-Gaussian corrections with the excursion set 
method, and we present our results for the halo mass function. 
We will see that, in the approximation in which the three-point 
correlator at different times is replaced by the corresponding 
cumulant, we recover eq. ( fTOb exactly, including the replace- 
ment H2i , except that now this replacement is not performed 
ad hoc to fit the A^-body simulations, but is the consequence 
of the diffusing barrier model of paper II, which also pre- 
di cts yja ~ 0.8 9, in remarkable agreement with the findings 
of iGrossi et ail J2009h . 8 We will see however that this result 
comes out in a rather unexpected way. In fact, the "local" 
term that, in excursion set theory, is supposed to give back 
the PS result multiplied by the appropriate factor of two, ac- 
tually vanishes, because for the three-point correlator (as well 
as for all odd-point correlators) the contribution of the image 
gaussian cancels the Press-Schechter contribution rather than 
adding up. The result ( fTOb comes entirely from non-trivial 
memory terms, that have no correspondence in the naive PS 
approach. 

We will then go beyond the approximation in which the 
three-point correlator at different times is replaced by the cor- 
responding cumulant, by computing explicitly the mass func- 
tion at next-to-leading order and at next-to-next-to-leading or- 
der in the small parameter <r 2 /S 2 . We wil then find further 
corrections, which depends on the derivative of the correlator, 
and which, with respect to the small parameter o 2 jS 2 , are of 
the same order as the subleading terms given in eq. ( flOb . 

Finally, in Section [4] we present our conclusions, summa- 
rizing the findings of this series of three papers. 

The focus of this paper is on the generalization of excur- 
sion set theory to non-Gaussian fluctuations. However, in 
appendix A we examine, with our path integral formalism, 
the generalization of naive PS theory to non-Gaussian fluctu- 
ations, and we will contrast it with the generalization of ex- 
cursion set theory. 

We have attempted to write this paper in a reasonably self- 
contained manner, but the reading of this paper will certainly 
be facilitated by a previous acquaintance with the first two 
papers of this series, in particular with paper I. 

passage of a barrier by trajectories performing a non-markovian stochastic 
motion (at least order by order in the non-markovian corrections). Of course, 
one should not forget that excursion set theory itself gives only an approxi- 
mate description of the physics involved. 

8 The parameter that we denote by a is the same as the parameter q of 
IGrossi et alj <2009t) . 



2. PATH INTEGRAL APPROACH TO STOCHASTIC 
PROBLEMS. NON-GAUSSIAN FLUCTUATIONS 

2.1. General formalism 

In this section we extend to non-Gaussian fluctuations the 
path integral approach that we developed in Section 3 of pa- 
per I for gaussian fluctuations. Our notation is as in paper I. 
In particular, we consider the density field 5 smoothed over a 
radius R with a tophat filter in coordinate space. We denote 
by 5 the variance of the smoothed density field and, as usual 
in excursion set theory, we consider S as a variable evolving 
stochastically with respect to the "pseudotime" 5, (see e.g. 
Sections 2 of paper I). The statistical properties of a random 
variable S(S) are specified by its connected correlators 

(6(S0...5(S P )) C , (13) 

where the subscript c stands for "connected". We will also 
use the notation 

(6 p (S)) c = [i p (S), (14) 

when all arguments 5i ,52, ■ ■ ■ are equal. The quantities /j. p (S) 
are also called the cumulants. As in paper I, we consider an 
ensemble of trajectories all starting at 5o = from an initial 
position (5(0) = So (we will typically choose So = but the com- 
putation can be performed in full generality) and we follow 
them for a "time" 5. We discretize the interval [0,5] in steps 
A5 = e, so = ke with k= l,...n, and 5„ = 5. A trajectory is 
then defined by the collection of values {Si,... , 5,,}, such that 
S(S k ) = S k . 

The probability density in the space of trajectories is 

W{So;d u ...,6 n ;S n )= (S D (S(Si )-Si)... 6 D (S(S„) - S„)) , 

(15) 

where Sd denotes the Dirac delta. As in paper I, our basic 
object will be 

n £ (5 ;4;5„)= / dS\... / dS n -iW(So;5 u ...,S n -uS n ;S n ). 



(16) 

The usefulness of n e is that it allows us to compute the first- 
crossing rate from first principles, without the need of postu- 
lating the existence of an absorbing barrier. In fact, the quan- 
tity 

dS n U e (S ;S n -S n ) (17) 



gives the probability that at "time" 5„ a trajectory always 
stayed in the region 6 <S C , for all times' smaller than 5„. The 
rate of change of this quantity is therefore equal to minus the 
rate at which trajectories cross for the first time the barrier, so 
the first-crossing rate is 

T(S ^ = -7^- I 1 d5„TlA5o;5 n ;S n ). (18) 

The halo mass function is then obtained from the first- 
crossing rate using eq. (fTJ together with (see eq. (33) of pa- 
per I) 

f(a) = 2a 2 F(a 2 ), (19) 
where 5 = a 1 . For comparison, it is also useful to introduce 



rips, e (c>o;'5«;5„) = / dSi...dS„-iW(5o;S u ...,S„-i,S„;S„). 

(20) 

So, Hp5, e (So;S n ;S„) is the probability density of arriving in S„ 
at time 5„, starting from So at time 5o = 0, through any possi- 
ble trajectory, while H e (So;S n ;S n ) is the probability density of 



arriving in 8„ at time S, again starting from So at time Sq = 0, 
through trajectories that never exceeded 5 C . Observe that in 
both cases the final point S„ ranges over — oo < S„ < oo. In- 
serting eq. (T5[ into eq. ( f20b and carrying out the integrals 
over dSi . . . d5 n -\ we see that 

Tl vs ^5 ;5 n ;S n ) = {8 D {5{S n )-5 n )). (21) 

Therefore Ilps.e can depend only on the correlators (fT~3T > with 
all times equal to S„, i.e. on the cumulants ^ p {S n ). In contrast, 
n e (^o; 8 n ;S n ) is a much more complicated object, that depends 
on the multi-time correlators given in eq. (TT~3T >. 

Furthermore, we see that Ilps.e is actually independent of e, 
since the integration over the intermediate positions has been 
carried out explicitly, and the result depend only on S„ and 
S n . Thus, we will write IIps ie simply as ITps. In contrast, IT 
depends on e, and we keep this e dependence explicit. We are 
finally interested in its continuum limit, II e= o, and we have 
already seen in paper I that taking the limit e — * of II e is 
non-trivial. So, despite their formal similarity, II £ and lips are 
two very different objects. The distribution function lips has 
a trivial continuum limit, and depend only on the cumulants, 
while n e depends on the full correlation functions ( Tl3l . and 
its continuum limit is non-trivial. All the complexity enters 
in Il e through the presence of a boundary in the integration 
domain, since the variables 8j are integrated only up to 8 C . 

The use of lips generalizes to non-Gaussian fluctuations 
the original PS theory, since we are integrating over all tra- 
jectories, including trajectories that perform multiple up- and 
down-crossings of the critical value 8 C , and therefore suffers 
of the same clou d-in-cloud problem of the original PS theory . 
In the literature (IChiu et all 1 19971: Robinson & Baker 1 120001: 
iMatarrese et alj|2000t lLo Verde et alj|2008l) this density func- 
tional has then been used together with the ad hoc prescription 
that we must multiply the mass function derived from it by a 
"fudge factor" that ensures that the total mass of the universe 
ends up in virialized objects. For gaussian fluctuations this is 
the well-known factor of two of Press and Schechter, while for 
non-Gaussian theories it is different, although typically close 
to two. 

In contrast, IT e generalizes to non-Gaussian fluctuations the 
approach of the excursion set method, where the "cloud-in- 
cloud" problem is cured focusing on the first-passage time of 
the trajectory, and no ad hoc multiplicative factor is required. 
So, II e is the correct quantity to compute. From the compar- 
ison of II e and ITps performed above, we understand that the 
difference between the two is not just a matter of an overall 
normalization factor. As we have seen above, in lips all the in- 
formation contained in the correlators at different "times" get 
lost, since it depends only on the cumulants. The correlators 
at different time contain, however, important physical infor- 
mation. Recalling that the role of "time" is actually played 
by S(R), the correlators at different time are actually correla- 
tors between density fields at different smoothing scales R\, 
R2, etc., and therefore carry the information on the depen- 
dence of halo formation on the environment and on the past 
history. These informations are intrinsically non-markovian, 
which is the reason why IT e is much more difficult to compute. 
However, these correlations are physically very important, es- 
pecially when we study the non-Gaussianities, and are com- 
pletely lost in the extension of PS theory based on lips. For 
this reason, our real interest is in computing the distribution 
function IT e , while lips will only be considered as a bench- 
mark against which we can compare the results provided by 
IT. 



The first problem that we address is how to express 
IIps(<5o;(5;S) and H e (5o;8;S), in terms of the correlators of the 
theory. Using the integral representation of the Dirac delta 



f°°dX _ iXx 



(22) 



we write eq. ( fT5l ) as 



W(S ;S h . . . , S n ;S n ) = [°° ^ . . . ^ AA (e^ 
J-oo 2tt 2tt 

(23) 

We must therefore compute 

e 2 = w,) ). (24) 

This is a well-known object both in quantum field theory and 
in statistical mechanics, since it is the generatin g functional 
of the connected Green's functions, see e.g. IStratonovichl 
(11967b . To a field theorist this is even more clear if we de- 
fine the "current" J from -iX = eJ, and we use a continuous 
notation, so that 



e z = c 



if dSJ(S)6(S)\ 



(25) 



Therefore 



p=2 <i = l i,,= l 

1 {-if 

(5jSj) c + XiXjXu (Si5jS k ) c 

(-0 4 

+— jj-A, AyAjfcA/ (Si5j5 k Si) c + . . . , 



(26) 



where Sj = S(Si) and the sum over ... is understood. This 

gives 



W(S ;S l ,...,S„;S n )= / VX 



(27) 



exp <; i \&i + £ -Jj- £ ■ ■ ■ £ K ■ ■ ■ K (Su ■ ■ ■ 5i„ > 

1=1 p=2 1 ' h = \ i„=l 



where 



so 



VX: 



2lT 27T 



(28) 



IIps(<5o;<5«;S„) = / d8\...d& n -\ VX 



(29) 



ex p { ' £ A a + £ -rr £ ■ ■ ■ £ A/ > ■ ■ ■ A/ " (*««••■ s h 

i=l p=2 i\=\ i p =\ 



and 



U e (So;S n ;Sn) = / d8\...d8„-i I VX 



(30) 



exp <; iJ2 Xt8i+J2 -TT £ • ■ ■ £ A '. ■ ■ ■ A v ( s u ■ ■ ■ ^ 

/1 = 1 p=2 i=l i p =\ 
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2.2. Perturbation over the markovian case 

As it was found in the classical paper bv lBond et al.l d!991l) 
when the density S(R) is smoothed with a sharp filter in mo- 
mentum space it satisfies the equation 

d6(S) 



dS 



(31) 



where here S = cr 2 {R) is the variance of the linear density field 
smoothed on the scale R and computed with a sharp filter in 
momentum space, while 7](S) satisfies 

(i0iM5d) = 8(S 1 -S2). (32) 

Equations OTb and (l32t are formally the same as a Langevin 
equation with a Dirac-delta noise r/(S). In this case, as dis- 
cussed in paper I, 

(5(S i )6(Sj)) c = mm(S i ,Sj), (33) 

and for gaussian fluctuations, where all n-point connected cor- 
relators with n > 3 vanish, the probability density W can be 
computed explicitly, 



w gm (£ ; <$i,... A;S„) = 



l 



(27re)"/ 2 



exp ■ 



1 n ~ l 

i=0 



3d 2 



(34) 

where the superscript "gm" (gaussian-markovian) reminds 
that this value of W is computed for gaussian fluctuations, 
whose dynamics with respect to the smoothing scale is marko- 
vian. Using this result, in paper I we have shown that, 
in the continuum limit, the distribution function n e= o(<5; 5), 
computed with a sharp filter in momentum space, satis- 
fies a Fokker-Planck equation with the boundary condition 
H e=0 (8 c ,S) = 0, and we have therefore recovered, from our 
path integral approach, the standard result of excursion set 
theory, 



Ui: (S ;S;S)=- 



1 



-(«5-(5o) 2 /(2S) _ -(2S-8 -5) 2 /(2S) 



(35) 



For a tophat filter in coordinate space, we have found in pa- 
per I that eq. d33l is replaced by 



(5(Si)8(Sj)) c = mm(S h Sj) + A(S u Sj) , 



(36) 



where S is now the variance of the linear density field com- 
puted with tophat filter in coordinate space. We found that 
(for the ACDM model used in paper I) A(Si,Sj) is very well 
approximated by the simple analytic expression 



A(Si,Sj): 



Si(Sj-Si) 



(37) 



where Sj < Sj (the value for Si > Sj is obtained by symmetry, 
since A{S,,Sj) = A(Sj,S,)), and k(R) is given in eq. @. The 
term min(5', ,5' y ) in eq. d36l l would be obtained if the dynam- 
ics where governed by the Langevin equation eq. OTb . written 
with respect to the variance S computed with the tophat filter 
in coordinate space, and with a Dirac delta noise, and there- 
fore describes the markovian part of the dynamics. The term 
A(S,,Sj) = A,j is a correction that reflects the fact that, when 
one uses a tophat filter in coordinate space, the underlying dy- 
namics is non-markovian. Observe that the full two-point cor- 
relator fl36l l cannot be obtained from an underlying Langevin 
equation and, as a consequence, the probability distribution 
H e (So',8 n ;S„) does not satisfy any local generalization of the 
Fokker-Planck equation, see the discussion below eq. (83) of 



paper I. However, the formalism developed in paper I allowed 
us to compute H e (8o;5 n ;S„) directly from its path integral rep- 
resentation, 



6c 



dSi . . .d5„-i / T>\ 



U e (S ;S„;S„) ■ 



x exp U\ i 8 i --[mm(S i ,S j ) + A(S i ,Sj)]\ i \jj , (38) 

by expanding perturbatively in A(Si,Sj). The zeroth-order 
term simply gives eq. d3~5l . i.e. the standard excursion set 
result, with the variance of the filter that we are using. The 
first correction is given by 



U.f\S ;S n ;S n ): 



1 

dS\ . . . dd n -i - ^ Aijdidj 



(39) 



x f PA exp < i ^ A; Sj - — ^ min(S,- , Sj) A,- Xj 



- J2 A U "dSi... d8 n . { didjWXSo; S u ...,8 n ;S n ), 



where we used the notation 9, = d/d8j anf the identity Xe lXt = 
—id x e . This quantity has been computed explicitly in Sec- 
tion 5.3 of paper I, and the corresponding result for the halo 
mass function is given by eq. (fj). In this paper we will per- 
form a similar computation for the correction induced by the 
three-point function. 

3. EXTENSION OF EXCURSION SET THEORY TO 
NON-GAUSSIAN FLUCTUATIONS 

If in eq. (f30b we only retain the three-point correlator, and 
we use the tophat filter in coordinate space, we have 

U t (S ;S„;S„) = / 1 dSi. ..dSn-i / VX (40) 



x exp jj'A,^ - i [min(Si,Sj) + A ij ]X i X j + (<5,£ s 6 k ) A/Ay A* 

Expanding to first order, Ay and (SjSjSk) do not mix, so we 
must compute 

i " r Sc 

Uf\S ; 5„;S„) = --J2 ( W*) / d ^ • • ■ d ^-\ didjd k W sm , 
i.j.k=\ J -°° 

(41) 

where the superscript (3) in n^ 3) refers to the fact that this is 
the contribution linear in the three-point correlator. In prin- 
ciple the expression given in eq. ( |4TI ) can be computed using 
the formalism that we developed in paper I. In the continuum 
limit the triple sum over i,j,k in eq. (HTb becomes a triple in- 
tegral over intermediate time variables dSi,dSj,dSk, each one 
integrated from zero to S„, so the full result is given by a triple 
time integral involving (S(Si)S(Sj)S(Sk)), which is not very il- 
luminating. 

Fortunately, such a full computation is not necessary ei- 
ther. Remember in fact that the non-Gaussianities are partic- 
ularly interesting at large masses. Large masses correspond 
to small values of the variance S = a 2 (M). Each of the inte- 
grals over dSi,dSj,dSit must therefore be performed over an 
interval [0,S„] that shrinks to zero as S„ — > 0. In this limit 
it is not necessary to take into account the exact functional 
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form of (5(Sj)S(Sj)8(Sii)). Rather, to lowest order we can re- 
place it simply by (S 3 (S„)). More generally, we can expand 
the three-point correlator in a triple Taylor series around the 
point S, = Sj = Sk = S„. We introduce the notation 



Gf^iSn): 



d p d q d 
dSfdS 1 j Ts k 



-(5(SdS(Sj)6(S k )) 



(42) 



Then 



(SiSiMSjMSk)) = 

(-Y)P+q+r 



(43) 



E 

P,1,r=0 



p\q\r\ 



(S„ - Sd p (S„ - SjW n - S k ) r G^ q ' r) (S n ) . 



We expect (and we will verify explicitly in the following) that 
terms with more and more derivatives give contributions to 
the function f(cr), defined in eq. (|T), that are subleading in 
the limit of small a, i.e. for a/S c <C 1. So, we expect that the 
leading contribution to the halo mass function will be given by 
the term in eq. (l43l with p = q = r = 0. At next-to-leading order 
we must also include the contribution of the terms in eq. d43l 
with p + q + r = 1, i.e. the three terms (p = l,q = Q,r = 0), 
(p = 0,q= l,r = 0) and (p = 0,q = Q,r = 1), at next-to-next-to- 
leading order we must include the contribution of the terms in 
eq. (l43l with p + q + r = 2, and so on. 

Observe that, in a general theory, the functions ,q ' (JS n ) 
with different values of (p,q,r) are all independent of each 
other; for instance, 



(S(Si)S 2 (S n )) 



(44) 



is in general not the same as 



dS 



(S\S)) 



(45) 



s=s„ 



so Gf'°(S n ) cannot be written as a derivative of G^' U ' U) (S„). 

The terms ' (S„) in eq. d43l must all be treated as in- 
dependent functions, that characterize the most general non- 
gaussian theory (except, of course, for the fact that Gf' q ' r) (S„) 
is symmetric under exchanges of p,q,r). However, for the 
purpose of organizing the expansion in leading term, sublead- 
ing terms, etc., we can reasonably expect that, for small S„ 



(0,0,0)/ 



-i(p,q,r) 
T 3 



(46) 



i.e. each derivative d/dSi, when evaluated in S, = S„, gives 
a factor of order 1/S n . This ordering will be assumed when 
we present our final result for the halo mass function below. 
However, our formalism allows us to compute each contribu- 
tion separately, so our results below can be easily generalized 
in order to cope with a different hierarchy between the various 



(S«). 

3.1. Leading term 
The leading term in n (3) is 



(47) 

where the superscript "L" stands for "leading". This expres- 
sion can be computed very easily by making use of a trick 



that we already introduced in paper I. Namely, we consider 
the derivative of nf m with respect to S c (which, when we use 
the notation n| m (<5o;^ n ;5„), is not written explicitly in the list 
of variable on which nf m depends, but of course enters as 
upper integration limit in eq. (fT6b). The first derivative with 
respect to S c can be written as (see eq. (B8) of paper I) 



d 



— n? m ((5o; S*\S n ) = J2 I dSi... d5n-i diW°» 



(48) 



since, when d/dS c acts on the upper integration limit of the in- 
tegral over ddj, it produces W(5\, . . . , 6, = S c , ■ ■ ■ , S„;S„), which 
is the same as the integral of djW with respect to dSj from 
Si = — oo to Si = S c . Similarly 



85} 



nf m (c5o; 6 n ;S„) = J2 dSi... dS n ^ didjW 



gm 



(49) 



see eqs. (B9) and (B10) of paper I. In the same way we find 
that 



d5l 



n-l p8 L . 

U^iSo; 5 n -S n ) =Yj Ml '-- d5n ~ l d ' d AW sm ■ 
iJ,k=i J -°° 



(50) 

The right-hand side of this identity is not yet equal to the 
quantity that appears in eq. d47| i, since there the sums run up 
to n while in eq. d50l > they only run up to n-l. However, 
what we need is not really n^(c>o;<5„;S„), but rather its in- 
tegral over dS n , which is the quantity that enters in eq. ( TT~8b . 
Then we consider 

[ 5C d5 n n? L \So;5„;S„) = -Us? 1 ) 

J-oo 6 

" pSc 

x / d5\ . . .d8 n -\d8 n didjdkW 
i,j,k=i J -°° 
and we can now use the identity 



gm 



(51) 



/ dS 1 ...dS n - 1 d8 n d i djd k W sa 



''dS*. 
''85} 



d8i...d6 n -idS n W sm 

dS n nf m (^ ;(5„;5„), 



(52) 



d6 n U^ h \5 ;S„;S n ) = -^-^ [ ' dS„ns m (So;S n ;S n ) . 

(53) 



6 86} „ 

From eq. d35l ). setting for simplicity 5q = 



nS,(So = 0;<5„A) = 



1 



y/2-KSn 



-5 2 J(2S„) _ -(2S c -5„f/(2S n ) 



(54) 

Inserting this into eq. ( f53l > we immediately find the result in 
the continuum limit, 



J<5„n ( e 3 «(0;<5„;5„): 



(Si) 



3V2^sl /2 



sj 

(55) 
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We now insert this result into eqs. (fT8l and dT9b and we ex- 
press the result in terms of the normalized skewness 



S 3 (a) =^{S 3 (S)). 



(56) 



Putting the contribution of n (3 L) together with the gaussian 
contribution, we find 



m= - 

1 7T 



1/2 



^£ p S 2 JQ.a 2 ) 



(57) 



,2 r 



1 + 



6<L 



2JI 
a 2 



1 +■ 



d\na V c 



Remarkably, t h is agr ees exactly with the result obtained by 
lLoVerde"et al. (2008), performing an Edgeworth expansion 
of the non-Gaussian generalization of Press-Schechter theory, 
see eq. ( fTOb . 

However, the fact that a naive non-Gaussian generalization 
of PS theory gives the same result that we have obtained from 
the non-Gaussian generalization of excursion set theory (at 
least to leading order for small a/6 c ; we will see below that 
the subleading term gets corrections) is somewhat acciden- 
tal, as can be realized as follows. In the sum over i,j,k of 
didjdu in eq. (l47l i. it is useful to separate the contribution with 
i = j = k = n from the rest. Recall that in PS theory the up- 
per integration limit for the variables dd\, . . . ,dS„-i is +oo 
rather than 6 C (which reflects the fact that in PS theory one 
looks at the probability that, at a given smoothing radius, the 
smoothed density is above threshold, regardless of whether 
it was already above threshold for some larger smoothing ra- 
dius). If in eq. (|47T > we replaced the upper integration limit S c 
with +oo, a derivative d, with i < n would integrate by parts 
to zero. The terms where at least one of the indices i,j or 
k is strictly smaller than n therefore have no counterpart in 
PS theory. The term where all indices are equal to n, 

in contrast, are local terms, which depends only on the cu- 
mulants rather than on the correlators at different points, and 
that can have a correspondence with PS theory. In the gaus- 
sian case, it is just such a local term that gives back the PS 
result, together with the factor of two that in PS theory was 
added by hand. Formally, this comes from the fact that, in the 
gaussian case, the excursion set probability distribution is the 
difference between the original PS gaussian and an "image" 
gaussian, and these two terms give contributions that add up 
when computing the first crossing rate. 

In the case of the three-point correlator the situation is how- 
ever different. Denoting by n (3 ' La) the contribution to n (3 L ^ 
obtained by setting i = j = k = n in eq. d47l >. we have 



nfi a) (0;<5„;S„): 



i<^)9„ 3 nfr (0;5„;5„), 



(58) 



and therefore 

c 5c 



dS n n ( e 3 j La) (0; 5 n ;S n ) = - - (<5 3 ) [d 2 U^ (0; 5 n ;S n )] 



S„=8 C 



= 0. 

(59) 

This result is in a sense surprising. Since PS theory gives 
a wrong normalization factor, missing a factor of two in the 
gaussian case, and a factor close to two in the non-Gaussian 
case, what is done in the literature when one uses PS theory 
is to take the PS result and multiply it by hand by a factor 
of two (or, for non-Gaussian fluctuations, close to two), as- 
suming that this would come out from a proper treatment of 
the cloud-in-cloud problem, i.e. from excursion set theory. 



We see however that this is not at all the case. In excursion 
set theory II^™ is a difference of two gaussians, see eq. ([35l l, 
so all its derivative with respect to 5„ of odd order, evaluated 
in S„ = S c are twice as large as for a single gaussian, but the 
function itself, as well as all its derivative with respect to S„ of 
even order, evaluated in 5„ = S c , are zero, i.e. the contribution 
from the second gaussian cancels the first contribution, rather 
than adding up. Since in eq. d59l appears the second deriva- 
tive of II^™ in 6„ = S c , this term vanishes. We therefore see 
that the logic behind the use of PS theory for non-Gaussian 
fluctuations, namely (1): compute with a naive extension of 
PS theory to non-Gaussian fluctuations and (2): multiply the 
result by hand by a "fudge factor" ~ 2, assuming that it would 
come out from a solution of the cloud-in-cloud problem, is not 
justified. For the contribution linear in the three-point corre- 
lator (J 3 ), this "fudge factor" is actually zero, and the result 
comes entirely from terms with at least one derivative <9,- with 
i < n, which have no counterpart in a non-Gaussian extension 
of PS theory. Above we have computed the excursion set the- 
ory result performing at once the sum over i,j,k, using the 
trick given in eq. d52b . In appendix|B]we compute separately 
the terms in the sum over i, j,k with one or more indices equal 
to n, and we c heck th at the y give back eq. ( f55l >. 

In Sections [3.21 and [331 we will compute the corrections to 
eq. d57b to next-to-leading and to next-to-next-to-leading or- 
der. We also need to take into account that the barrier must 
be treated as diffusing, see paper II, and we must include the 
corrections due to the tophat filter in coordinate space. This 
will be done in Section [3~4l 

Before leaving this section we observe that, in the approxi- 
mation in which the correlators are replaced by the cumulants, 
the effects of the higher-order correlators can also be com- 
puted very simply. For instance, the effect of the four-point 
function (5*) is obtained using 
r$ c 

dSi... dd„-i d8 n didjd k d,W sm 



E 



d6 n U^iSo-J^S). 



(60) 



3.2. The next-to-leading term 
Using eqs. d4TT > and d43l . at next-to-leading order we get 

1 

2 



d6 n Il^ NL \do;S n ;S n ) = - G^'°'°\s n ) (61) 



x^(S„-Si)E / dSi...dS n . 1 d5„d i d j d k W sm , 
i=i j,k=i 00 

where the superscript "NL" in n ( e 3 ' NL) stands for next-to- 
leading, and we used the fact that the three terms (p = 1 , q = 
0,r = 0), (p = 0,q= l,r=0) and (p = 0,q = 0,r = 1) give the 
same contribution. We now use the same trick as before to 
eliminate J2"jk=i ^jdk in favor of d 2 /dS 2 , 

^ d8 n Uf^\So;S n ;S n )= l -Gf' '°\s n ) 

" d 2 f Sc 

X Y^(S n -Sd m / dSi . ..dSn-idSnW*" . 

The remaining path integral can be computed using the tech- 
nique developed in paper I, namely we write 



d6i...d5*-id6 n d i W eai 



(62) 
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= / d6i . . .dS„-idS„W(5o;S u . . . ,Si = 8 C ,. . . ,8 n ;S„) , 

J -co 

and we use 

W sm (5 ; S u ■ ■ ■ , Si-i , S c , 5 M S n ;S n ) (63) 
= W° m (6 ; ft, ... , 8m , 5 c ;Si)W^(5 c -5 M ,...,6 n ;S n -Sd, 

so 

rS e rS c 

ddi . . . dSi-i / dSi+i . . . dS„-idS„ 



x W gm (5 ( >; 8i,..., 5m , S^SiW^iSc, 6 M ,...,6„;Sn-Sd 
= m m (S ; ft; S t ) [ ' dS n IB m (ft; S„;S„- S,) . (64) 



Recalling from paper I that 

nf^ ;^;S) = V^^7^ 



<S c -8or/(2S) 



+ 0(e) (65) 



and 



1 S c -S„ 



7T S 3 / 2 



-(8-5„yl(2S) 



+ 0(e), (66) 



we see that the factors y/e in nf m (ft;ft;S) and in 
nf m (ft;ft;S) combine with ^ ( to produce an integral over 
dSi, and 



dd n n^ L \6o;6„-,S n )=^G ( t m (S n ) 
Air - 



dS, ■ 



o s 3 / 2 (s n -Siy/ 2 
a 2 



(67) 



dS 2 



f dS n {S t ft)exp{ 



The integral over d8 n is easily performed writing 
(ft-ft)exp^ -— — — V = (S n -Si)d„exp 



2(S n -Sd 



2(S n -Sd 



(68) 

so it just gives (S„-Sj). Carrying out the second derivative 
with respect to ft and the remaining elementary integral over 
dSi we get 

d6 n U^(S ;S„;S n ) = -±= ° 3 J n) e^^ . (69 ) 



-oo 

We now define 



U 3 (<t) : 



2tt S„ 
3G ( > m (S) 



(70) 



where as usual S = a 2 . When the ordering given in eq. d46*| > 
holds, U 3 (a) is of the same order as the normalized skewness 
83(a) given in eq. d56b . Computing the contribution to f(<j) 
from eq. (l69l and putting it together with eq. ( |57] | we finally 
find 



/(*) = 



2 \ 1/2 5, 



1 + T7-%G(0") 

6ft 



(71) 



where 



ft 4 



h NC (a)=^-S 3 (a)-^ 2S 3 (a) + U 3 (a) 



dS 3 
din a 



63(0-) + 6/3(0-)+ — !-■ 



oflnc dlna 



(72) 



We have ordered the terms in /zng(c) according to their im- 
portance in the limit of small cr/ft assuming, according to 
eq. 146b . that U 3 {a) is of the same order as 53(17). The lead- 
ing order is given by {5 C / cffSiia) and, as we have seen, it 
comes only from n (3 L) . The next-to-leading order in /zng(c) 
is given by the terms proportional to (ft/ a) I , and we see that 
it is affected by the terms with p + q + r = 1 in the expansion 
of eq. l43"T l. The terms in /z N g(ct) which are 0(1) with re- 
spect to the large parameter ft/er are next-to-next-to-leading 
order corrections and, if we wish to include them, we must for 
consistency include also the contribution from the terms with 
p+q + r=2 in the expansion of eq. 1431 . We compute them in 
the next subsection. 

Observe also that typically £3 depends very weakly on the 
smoothing scale R and hence on a. For instance, in / NL - 
theories it changes only by a factor ~ 3 as R is changed by a 
factor 100, from 0.1 Mpc/h to lOMpc/ft see lMatarrese et alj 
(2000). Therefore, even if parametric ally dS-^/d In a has the 
same power-law behavior as S3, its prefactor will typically be 
numerically small. 

3.3. The next-to-next-to-leading term 

Using eqs. d4TT > and 143b and keeping the terms with p+q+ 
r = 2 we find two kind of contributions. The first has (p = 2,q = 
r = 0), with a combinatorial factor of three and the second has 
(p = q = 1 , r = 0), again with a combinatorial factor of three. 
We denote the contribution to n (3) at next-to-next-to-leading 
(NNL) order by n (3 ' NNL) , and the two separate contribution 
with (p = 2, q = r = 0) and with (p = q = 1, r = 0) as n (3 NNLa) 
and n (3 NNLb) , respectively. Thus, 



"4 G3 



(2,0,0) 



(Sn) 



i=l 



dftnf< NNLa) (ft,;ft;ft,) 

(S n -Sd 2 J2 / d5i...d8 n - l d8 n d i d j d k W^, 

: J._l J -OO 



(73) 



7>1 ' 



and 



/ ' d5 n nf^\5M n ) = -\G^ fi \S n ) (74) 
J -00 



k=\ 



d8i...d5 n - l d5 n d i djd k W gn 



xJ2(Sn-Sd(S„-Sj)] 

The first term is straightforward to compute. We use again 
the trick of eliminating YTjk=\ ®i® k m f avor °f d 2 /dS 2 
we proceed just as in Section [3~2l The result is 

[ SC dS n nf^\5 ; S n ;S n ) = -l- Gf m (S n ) 



and 



2ttS 



l/2 e -Si/(2S„)_ 



ttS c Erfc 



/2Sn 



(75) 



where Erfc is the complementary error function. 

The computation of eq. ( |74T i is more complicated, but can 
be performed with the formalism that we have developed in 
paper I, see in particular appendix B of paper I. The factor 
SLi ®k is eliminated as usual in favor of d/dS c . We also ob- 
serve that, in eq. (|74l ), the terms with i = n or j = n do not 
contribute, because of the factor (S„-Si)(S„-Sj), and we sep- 
arate the sum over i,j into the term with i = j and twice the 
term with i < j. The first term is 

n-l 



i=l 



10 



'as. 



r s c 

■ dSi...dS n -idS„dfW gm , (76) 

J—oo 



and the second is 



n—2 n— 1 



i=i j=i+i 



d_ 

dS c 



d6 l ...dS n -idS n d i d j W sm . 



(77) 



As we discussed in detail in paper I, quantities such as the 
right-hand side of eq. d74l i are finite in the continuum limit 
e — > 0, as it is obvious physically, and as we checked explic- 
itly in solvable examples in paper I. However, when we split 
the sum over the indices z, j into two separate parts, such as 
those given in eqs. dToT l and ( 1771 ). these are separately diver- 
gent in the continuum limit, and the divergence cancels when 
we sum them up. It is therefore necessary to regularize them 
carefully, and separate them into a divergent part and the finite 
part. Since we know that the divergent terms must cancel, we 
can simply extract from each term the finite part, disregarding 
the divergences. This is the finite part prescription discussed 
and tested in detail in paper I. We will denote by TV this pro- 
cedure of extracting the finite part from terms such as d76b and 
d77]). 

The computation of the finite part of ( fToT l is basically iden- 
tical to the one that we already performed in appendix B of 
paper I, see in particular eqs. (B13)-(B15) and (B29) there, 
and the result is that this term diverges as 1/ ^/e with no finite 
part, so 

-1 a P S C 



_d_ 



dSi... dS n -xdS n dfW m = . (78) 



The computation (|77| i is also completely analogous to the 
computation of the "memory-of-memory" term performed in 
appendix B of paper I, see in particular eqs. (B17)-(B28) and 
(B30) there, except that we now have a factor (S„-Si)(S n —S /) 
in the integrals over dS, and dSj. We can then repeat basically 
the same steps as detailed in appendix B of paper I, and we 
find 



7rv27T 



TV 



dSi 



(S„-Sd 

c 3/2 



Si 



dSj 



(s„-sjy/ 2 



exp 



-s;/2s t 



(79) 



, Si (Sj-Si)'/ 2 -~ r { 2(Sj-Sd, 

where a = ^/ae and a is a numerical constant which appears 
when the sum over j is replaced by an integral over dSj, see 
eqs. (B20)-(B24) of paper I. The integral over dSj can be com- 
puted writing t„ = S„ - Si, tj = Sj - St, and using the identity 

(t n -t^ 2 < ~ 2 



-JjT 

t V2 



V/27T 

2a 



2 l/2 e - fl V(2r„)_ fl ^ Erfc 



I2t n 



(80) 



which is proved in the same way as eqs. (115) and (116) of pa- 
per I. In this equation a = y/ae goes to zero in the continuum 
limit. In the limit a — > the above result displays a term diver- 
gent as 1 j a, i.e. as 1 / ^/e, which must cancel the divergence 
coming from d76l >. plus a term which is finite as a — > 0, which 
can be extracted from eq. d80l ) recalling that Erfc(0) = 1, so 



TV 



dSj 



1/2 



(S n -Sj) 



exp 



2(Sj-Sd 



(81) 



Computing the remaining integral over dSi, which is finite and 
elementary, we find 

f l dS n Uf^ b \S ;S„;S n ) = -I G^ m (S„) 



So 



r^s\! 2 e<l^ -TrcLErfc 



Putting together this result and eq. (F75T > we end up with 



(82) 



^„n^ NNL '(<5 ;^;5„) 



2tt \2 



(1,1,0), 



1/2 -S 2 C /(2S„)_ 



2irS z e 



We now introduce the function 



TrS r Erfc 



V 3 (£r) : 



(83) 



(84) 



According to eq. d46| i, Vj(a) is parametrically of the same 
order as 6>3(er) and U}{(t), as a — > 0. We can now compute 
the contribution of this term to the function /ing(c) using 
eqs. ( fT8l . ( fT9] l and dTTI l. Retaining only the terms that con- 
tribute up to 0(1) in S c / u, we find that the function /ing(c) is 
modified to 



c4 si / jc 

/*ngM=-t«S 3 (ct)— I [2S 3 (a) + Lk(a)-— — 
o a L \ dm a 



d\na d\na 



S 2 



(85) 



This is the complete result for the halo mass function, up to 
NNL order in the small parameter a 2 /5 2 . 

3.4. The effects of the diffusing barrier and of the filter 

Until now we have worked with a barrier with a fixed height 
S c and we neglected the corrections due to the filter. We now 
include the modifications due to the fact that the height of the 
barrier diffuses stochastically, as discussed in paper II, and 
also the corrections due to the filter. 

To compute the non-Gaussian term proportional to the 
three-point correlator with the diffusing barrier we recall, 
from paper II, that the first-passage time problem of a particle 
obeying a diffusion equation with diffusion coefficient D = 1, 
in the presence of a barrier that moves stochastically with dif- 
fusion coefficient Db, can be mapped into the first-passage 
time problem of a particle with effective diffusion coefficient 
(I+Db), and fixed barrier. This can be reabsorbed into a 
rescaling of the "time" variable S — > (1 +Db)S = S/a, and 
therefore a — > cr/y/a. At the same time the three-point cor- 
relator must be rescaled according to (<^) — > a~ 3 / 2 (<^) since, 
dimensionally, (5^) is the same as 5 3//2 (if we perform di- 
mensional analysis as discussed below eq. (A10) of paper I), 
which means that 53 — » a I//2 6>3, and similarly for the functions 
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Ut, and V3. 9 Then eqs. ( |7TT i and (l85l l become 



m= - 

7T 



2\ 1/2 fl '/M 



£ e -fl5 c i /(2<r i ) 



1 + 



6a ^ 



■%g(o-) 



(86) 



where 



,2x4 



a5l 



din a 



S 3 (cr) + U 3 (cr) + V 3 (cr) + 4r^~ + ^ ^ 



ii In cr diner 



+0 



(87) 



We see that the terms depending on the skewness 6>3(cr) and its 
derivative coincide with those given in eq. ( [Tol l, if we identify 
Jeff with a l l 2 8 c . Observe, from eq. (Q, that our prediction 
a 1 / 2 ~ 0.89 is in r emarkable ag r eemen t with the value a 1 / 2 ~ 
0.86 proposed by iGrossi et all (120091) from the fit to the N- 
body simulations (see however footnote^. 

We have therefore derived, from a fi rst principle computa- 
tion, e q. dTOb. which was proposed in lLo Verde et al.l (2008) 
and in Grossi et al. (2009) using a mixture of heuristic the- 
oretical arguments (the use of a non-Gaussian extension of 
PS theory, rather than of the excursion set theory) and a cal- 
ibration of parameters from the fit to the data of the A^-body 
simulations (the replacement S c — > 0.86<5 e ), and we have im- 
proved it including th e effect of the f unctio ns U^{a) and Vi(a), 
which are absent in lLoVerde e t al. (2008) and cannot be ob- 
tained from any naive extension of PS theory, which from the 
beginning contains only the cumulants, rather than the full 
correlation functions at different smoothing radii. 

The term in eq. (f87t which is dominant for small a is the 
same as that of both eqs. ([8]) and dTOb . a nd appears to fit well 
the data of the A/-body simulations (Gro ssi et alj20 09). Given 
the size of the error bars of the n on-Gaussian jV-body simu- 
lations (see e.g. Fig. 6 and 7 of IGrossi et al.l (|2009)), it is 
probably difficult for the moment to test the subleading terms 
in eq. dSTT i, and in particular to see the effect of the functions 
U 3 (a) and V 3 (cr). 

As a final ingredient, we must add the effect of the 
tophat filter function in coordinate space. When the non- 
gaussianities are not present, these are given by eq. ©■ More 
generally, even the non-Gaussian corrections must be com- 
puted using the propagator [min(Si,Sj) + Ay] in eq. d40b . so 
we will apply the same correction factor found for the gaus- 
sian part also to the non-Gaussian term, and we end up with 



f(a) = (l-R) - 

7T 



1/2 a^Sr 



-aS 2 /(2a 2 ) 



a 



1 + 



6^/2 A 



■h NC (a) 



+ ^^r(o, o , 

sjl^ a V 2cr 2 



(88) 



with /7ng(c) still given by eq. d87l >. More generally, also the 
term proportional to the incomplete Gamma function could 
get non-Gaussian corrections, which in principle can be com- 
puted evaluating perturbatively a "mixed" term proportional 

9 In principle we should also shift the argument a of S3, U3 and V3. How- 
ever, 1S3 depends very weakly on the smoothing scale R and hence on a. For 
instance, in /nl -theories it changes only by a fact or ~ 3 as R is c hange d by 
a factor 100, from 0.1Mpc//i to 10Mpc//i, see Matar rese et alj <2000l) . In 
most situation, we can then neglect the rescaling of the argument of 1S3 , and 
we expect that the same holds for U3 and V3 . 



to 



Aij(S k Si8 m )didjdkdid m 



(89) 



in eq. (|40l . However we saw in paper I that in the large mass 
limit, where the non-Gaussianities are important, the term 
proportional to the incomplete Gamma function is sublead- 
ing, so we will neglect the non-Gaussian corrections to this 
subleading term. 10 

The relative weight of the correction due to the filter pro- 
portional to the incomplete Gamma function, and of the non- 
Gaussian corrections depends on the value of a and, of course, 
on the value of the three-point correlator, i.e. of S3. In 
/nl theory S3 increases very weakly with the mass, i.e. as 
a — > 0. In the low-er (i.e. large mass) limit we can use 
the asymptotic expansion of the incomplete Gamma func- 
tion for large z, T(0,z) ~ z~ l e~ z , and we see that, asymptot- 
ically, the term in the second line of eq. (l88l depends on a as 
crexp{-a<5 2 /(2cr 2 )}, and therefore is small compared to both 
the leading and next-to-leading term in the non-Gaussian cor- 
rections, which overall behaves as <j~ 3 St,((t) exp{-a<5 2 / (2c 2 )} 
and cr~ 1 iS 3 (cr)exp{-c!(5 2 /(2cT 2 )}, respectively. When this 
asymptotic behavior sets in depends, of course, on the nu- 
merical value of S3 so, in /nl -theory, on the value of the 
/nl parameter. In any case, given a measure of /(cr), either 
from galaxy surveys or from iV-body simulations with non- 
Gaussian initial conditions, the prediction (l88l allows us to 
disentangle the effects due to the filter from the physically in- 
teresting effects due to primordial non-Gaussianities. 

4. CONCLUSIONS 

To conclude this series of three papers, we summarize the 
main results that we obtained and, at the price of some rep- 
etition, we collect here the most important formulas that are 
scattered in the text. Our aim was to compute the halo mass 
function, i.e. the number density n(M)dM of dark matter ha- 
los with mass between M and M+dM, both for gaussian and 
non-Gaussian primordial density fluctuations. This can be 
written as 



dM 



M 2 dlnM 



and the issue is to compute the function f(u). Our final result 
can be written as 



f(a) = (l-R) 



2\ 1/2 « 1 / 2 <5, 



c_ -a8 2 J{2a-) 



a 



1 + 



6a i m t 



-/ing(ct) 



id;. 



2o l 



'lit cr 



(91) 



where T(0,z) is the incomplete Gamma function. Three dis- 
tinct physical effects are taken into account in this result. 

One is the fact that we have treated the threshold for gravita- 
tional collapse as a stochastic variable that fluctuates around 
an average value, which is S c ~ 1.686 for the spherical col- 
lapse model, and is a rising function of a for the ellipsoidal 
collapse model. As discussed in paper II, this is a way of 
taking into account, at least at an effective level, part of the 
complexity of a realistic process of halo formation, which is 
missed in the simple spherical or ellipsoidal collapse model. 

10 Furthermore, one must be aware of the fact that the term propor- 
tional to /inq(o") might in general receive corrections from the tophat filter 
that have not exactly the same form as that of the gaussian term. Again, 
in principle these can be obtained by computing the term proportional to 
{&k$i$m) did jdkdid,,, in the expansion of the path integral. 
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Furthermore, the stochasticity of the barrier reflects uncertain- 
ties in the operative definition of what is a dark matter halo. 
The inclusion of a diffusing barrier gives rise to the constant a 
in the above result. This constant enters also in the exponen- 
tial, thereby modifying dramatically the behavior predicted by 
PS theory. Our prediction is a ~ 0.80, i.e. y/a ~ 0.89, which 
gives a remarkable agr eement with t he da ta from A^-body sim- 
ulations. For instance TGrossi et ail d2009l) . from the fit to the 
Af-body simulation, find y/a ~ 0.86. 

A second effect included in eq. d9TT > is that we have prop- 
erly accounted for the fact that the comparison with the data, 
whether observational or from A^-body simulations, requires 
the use of a tophat fi lter function in c oordinate space. In 
the classical paper of lBondet all (1 19911) . using a tophat fil- 
ter in momentum space, the computation of f(a) was reduced 
to a first-passage time problem for a quantity that obeys a 
Langevin equation, and therefore the underlying dynamics is 
markovian. When one considers a different filter function, the 
dynamics becomes non-markovian and therefore the problem 
is much more complicated. Basically, this is the issue that for 
a long time blocked further analytical progress on this prob- 
lem. In paper I of this series we have developed a formalism 
in which the problem is formulated in terms of a path integral 
with boundaries, and non-markovian corrections can be com- 
puted perturbatively. In eq. d9Tb this effect enters through the 
constant R, defined as k = ok with n given by eq. (0|. 

The third effect, which was the subject of the present paper, 
is the inclusion of the non-Gaussianities. These are contained 
in the function /!ng(c). Using the path integral technique 
developed in paper I, we have computed it to leading, next- 
to-leading and next-to-next-to-leading order in the parameter 
o 2 /SI, which is small for large halo masses, where one can 
hope to see the effect of non-Gaussianities on the halo mass 
function. Our result is 

g a L \ dm a 



dlna diner 



SI 



(92) 



The functions S 3 (a), U 3 (a) and V 3 {a) are defined in terms 
of the three-point correlator of the smoothed density field 
(S(S\)S(S2)S(Si)) and of its derivatives, as follows, 



S 3 = ±(6 3 (S)), 



U 3 = - 



dS 



■{5{S,)5\S)) 



(93) 
(94) 



do i 



Si=S 



s,=s 



+12 



d d 
dS\ dSz 



and we prefer to write them as functions of a = y/S, 

Our result has passed to a good accuracy various compar- 
isons with numerical results. First of all, one can study nu- 
merically what happens in the excursion set theory, with fixed 
(rather than diffusing) barrier and tophat filter in coordinate 
space, by performing a Monte Carlo realization of the first- 
crossing distribution of excursion set theory, obtained by inte- 
grating numerically a Langevin equation with a colored noise. 
This w as recent ly performed in d etail in Ro bertson et alJ 
(2008) (see also iBond et al.1 (119911) ). In this limit our ana- 
lytical result is obtained from eq. ( |9TT > setting a = 1 (since the 
barrier is taken as fixed in the Monte Carlo simulation) and 
^ng(c) = 0, i.e. we are testing the effect of k. Comparing 
our results in paper I with Fig. 4 of Robertson et al. (2008) 
we find very good agreement. This is a first useful test of our 
technique. 

Using eq. d9TT l with a ~ 0.80 (obtained by reading the dif- 
fusion coefficient of the barrier Db from Af -body simulations, 
and using our prediction a = 1 /(l +Db)) and with k given in 
eq. dU, and setting /ing(o') = 0, we can compare our result with 
the mass function found in A^-body simulations with gaussian 
initial conditions. The comparison is shown in Figs. 6 and 7 
of paper II. For all values of er" 1 > 0.3 the discrepancy be- 
tween our analytic result and the Tinker et al. fit to the same 
A^-body simulation is smaller than 20%, and for a~ l > 1 it is 
smaller than 10%. Considering that our result comes from an 
analytic model of halo formation with no tunable parameter 
(the parameter a is fixed once Dg is given, and we do not have 
the right to tune it), while the Tinker et al. fitting formula is 
simply a fit to the data with four free parameters, we think that 
this result is quite encouraging. The numerical accuracy is ac- 
tually the best that one could have hoped for, considering for 
instance that we have neglected second-order non-markovian 
corrections. 

Finally, our prediction for the function /!ng(c) can be tested 
against A^-body simulations with non-Gaussian initial condi- 
tions. To leading o rder in the smal l a lim it, o ur result reduces 
to that proposed bv lLoVerde et alJ {2008) and Matarres e et alJ 
(2000) using non-Gaussian exte nsions of PS theory, and it has 
been found in lGrossi et al. ( 2009) that this formula reproduces 
very well the data, see in particular their Figs. 6 and 7. The 
size of the error bars is probably still too large for discrimi- 
nating between different forms of the subleading term. 

We thank Sabino Matarrese for useful discussions. The work 
of MM is supported by the Fond National Suisse. The work 
of AR is supported by the European Community's Research 
Training Networks under contract MRTN-CT-2006-035505. 
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(95) 
APPENDIX 



A. EXTENSION OF PRESS-SCHECHTER THEORY TO NON-GAUSSIAN FLUCTUATIONS 

As we repeatedly emphasized, the really interesting quantity for comparison with experimental data from galaxy surveys, 
and with A^-body simulations, is the distribution function II e , that generalizes excursion set theory to non-Gaussian fluctuations. 
The function lips defined in eq. ( 1201 ). where the integrations over the variables ddj run up to +oo rather than up to 6 C , not only 
suffers from the fact that it predicts that only a fraction of the total mass of the Universe finally ends up in virialized objects 
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(the infamous factor of two that Press-Schechter where forced to introduce by hand) but also misses all the subtle correlations 
between different scales which are just one of the characteristic features of non-Gaussianities. For this reason, in the body of 
this paper we concentrated on the computation of II e . Still, it is interesting to see how our path integral formalism reproduces 
PS theory and generalizes it to non-Gaussian theories. We discuss the issue in this appendix. In particular, we will see that, 
even in the non-Gaussian case, ITps satisfies a differential equation which is local in "time", the Kramers-Moyal equation, and 
which generalizes the Fokker-Planck equation. It is interesting to contrast this result with what happens for IT £ which instead, as 
discussed in paper I, does not satisfy any local diffusion-like equation. 

With our "microscopic" formalism based on the path integral, it is very easy to derive PS theory and to extend it to non- 
Gaussian fluctuations. Simply, in eq. (l29l l each integral over dSt, with 1 < i < n— 1, produces a factor 2tt8d{\), which allows 
us to perform trivially all the integrals over d\i with i < n. Denoting the residual variable A„ by A and setting for notational 
simplicity So = (the general result is recovered with 5 — > S — So), eq. (|29l l becomes 



U FS (S = 0;S;S) = 



— exp{/A(5 + ^2 —V P (S)} ■ 

p=2 P - 



When all ji p with p > 3 vanish, the integral gives a gaussian and we get back the standard PS result, 



n ps ((5o = 0;<5;S): 



1 



(2ttS) 1 / 2 



-6 2 /2S 



(Al) 



(A2) 



since, by definition (12(8) = S, where S is the variance computed with the filter function of our choice. Equation dAU generalizes 
PS theory to arbitrary non-Gaussian theories. 11 Observe that eqs. ( IA1 b and ( IA2b hold independently of the filter function used, 

and the fi p are the cumulants computed with the filter function in which one is interested. 

Equation ( IAU is a well-known re sult in the theory of stochastic processes (see e.g. I Risked d!984l) ). and it was applied to 
/NL-theory in Mat arrese et al.l d2000t) . Using this expression, the usual strategy in the literature is to compute J-(S) using 



d_ 
dT 



dxUps(S ;S;S), 



(A3) 



and to multiply by hand by a fudge factor ~ 2 to ensure the proper normalization. As we have shown in the discussion below 
eq. d59b . this multiplication by a fudge factor is not justified for non-Gaussianities. Still, let us discuss from the mathematical 
point of view the properties of the function lips, in order to contrast them with the excursion set theory distribution function II e . 

First of all, it is instructive to rederive the expression ( IA1 b for lips, with generic filter and generic non-Gaussian theory in an 
alternative way, using the technique developed in paper I for computing the effect of the correction Ay to the two-point function, 
see eqs. ( f38l > and (f39b . To compute IIps(<5o;<5;S') when the two-point correlator (5iSj) c is generic, rather than equal to min(5',', Sf), 
and in the presence of the higher-order correlators, we write 



Observe that eB„ 



(8i8j) c = min(S i) ^) + [(^) c -min(S i ,Sj)] = eA u + tB u . (A4) 
= fJaQS)— S=0. We then expand the exponential in eq. d29| i in powers of eB,j and of the higher-order correlators, 



n P s(<So;<5„;S„) = 



d8\ . . . d5 n -\ I T>\ 
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1 + 9 J2 eBi i did i ( s i s A)cdidjd k + . 



'■../= 1 



ij,k=i 



JvXe'^'^-i TTij^Mi^ (A5) 



where di = d/dSj. The derivatives 9, contribute only when the index i = n, otherwise we have a total derivative with respect to 
an integration variable, and the corresponding boundary terms at 5 = ±00 term vanish. Here it is crucial that one integrates up to 
+00. When we rather consider II e , instead of lips, the upper integration limit is S c and we remain with complicated and non-local 
boundary terms, compare e.g. with eq. (83) of paper I. For lips however this boundary term is absent and 



n P s(<5o;<5«;S„) = 



1 1 r°° 

1 " ^ (S 3 „)c dl + . . . J dS,... dSn., W sm (S ; S l ,...,5 i ,...,S„;S„) 
n°- gau ((5o;<5;S). 



(A6) 



" A word of caution is necessary when one considers eq. )AU with correlators fi p with p > 4. For instance, keeping only M3 an d /M> one i s faced with 
an integral that diverges, since /i4(5) = (<5 4 (5)} > 0. The correct statement is that lips (<5o <5; S) is given, order by order, by the expansion of eq. )AU in powers of 
/14. However, the expansion in powers of /14 is only an asymptotic series, which can be used to approximate the true result up to a finite order in /14, but diverges 
if we keep an infinite number of terms. If instead the highest cumulant that we include in eq. )AU is the integral converges because (— i) 6 M6 = — A*6 < 0, while 
the integral diverges again if the highest cumulant that we include in eq. )AU is fig, since (-i) 8 M8 = +M8 > 0> an d so on. Anyhow, the whole issue of the full 
resummation of the contributions of the /14 or higher-order correlators is physically irrelevant. These correlators are in general computed using phenomenological 
parametrization of the non-Gaussianities, such as /NL-theory, that are meant to be a useful description of the true non-Gaussianities only to leading, and at most 
next-to-leading order in /nl, so in general only the first few terms in the series makes sense physically. 
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Since the derivative d„ = d/dS„ does not act on the correlators (which are functions of S„, but not of 5„), the expansion in 
the square brackets can be exponentiated back, and we can write 

n PS (<5o; S;S) = e^lJ ^ u (8 ;S;S) , (A7) 

where (using now Sq generic) 

lf^(d ;6;S) = -J—e^V™, (A8) 

(2-Ksyi 1 

and the differential operator K^g is given by 

p=3 

To prove the equivalence of eqs. (IA7l i and (IAU we write eq. (lAll i as 

/°° dX f 1 1 
— exp|/A<5-- Ai2 (5)A 2 + W NG (A)| , (A10) 

with 

Wng(A) = ^ — p-/x P (5), (All) 

and we expand the exponential in powers of Wng(A). Using X p e' Xx = (—id x ) p e we see that W NC (X)e iXx = W(-id x )e iXx , and the 
same holds for any power of Wng(A), so 

exp{W NG (A)} e /A * = exp{W NC (-id x )}e iXx = e k ™e iXx . {All) 

Therefore eq. ( lAlOt becomes 

/°° dX ( 1 1 
— exp|/A<5-- M2 (5)A 2 | , (A13) 

which agrees with eq. (IA7I >. So, the distribution function lips that gives the extension of the PS formalism to non-Gaussian 
fluctuations can be written equivalently in the integral form dAU or in the differential form dA7| ), with K N q given by eq. JA9I ). 

It is interesting to observe that the function IIps(<5o = 0;S;S) obeys a local differential equation, both in the gaussian and in the 
non-Gaussian case. Consider first a gaussian theory with a generic filter function, so p, p (S) = for p > 3. In order to see exactly 
where enters the difference between integrating up to 8 C , as in II e and integrating up to +00, we start from the definition 

n PS {S -S„;S n )= I d5i...d5 n -x f VXe iXi5t -^ Si ^ x>x > ', (A14) 



where the sum over 2, j = 1, . . . ,n is understood, and we derive a differential equation satisfied by lips, by taking the derivative 
with respect to S„, 

l|f = ("5 ^Mr) ■ " d5 "~ l I VXX « X ^ x A iX &- 5^ A ^} • (A15) 

Again, using Aexp{/Ax} = -id x exp{iXx}, inside the integral we can replace Xk — > — idk and A; — > -/<9/. Since we integrate over 
d6i,...dS n -\, but not over t/<5„, \f k^n the term c^, when integrated over dJ^, is a total derivative and gives zero, because at 
the boundaries <5* = ±00 the integrand vanishes exponentially, and the only contribution comes from k = n. Similarly, also di 
contributes only when I = n. This is the step that does not go through for II e , when the integration is only up to 8 C , and a 
complicated boundary term arises, see eq. (83) of paper I. Therefore, since (5 2 ) = S„, we get a Fokker-Planck equation 

dIl PS 19 2 n PS 

^ = 2^' (A16) 
whose solution, on the line -00 < S < 00, is indeed given by eq. (IA21 >. Equation (IA16b can be generalized to the non-Gaussian 
case using the integral form of the solution (lAlb and taking the time derivative, 

^-t^i r^ HXY J iU ,±^ S) \ , (A17) 

where £i p = dp p /dS. Inside the integral we can replace (iX) p e by d P e' Xx , so 

p=2 r 

This equation is call ed the Kramers-Moyal (KM) equa tion or "the stochastic equation", and is well known in the theory of 
stochastic processes dStratonovichl(ll967l) . lRiskenl d!984l) ). 
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In conclusion we have seen that, independently of choice of filter function, n ps satisfies a local differential equation both in 
the gaussian and in the non-Gaussian c ase. In the gaussian case it satisfies the FP equation dA161 >. while in the non-Gaussian case 
it satisfies the Rramer-Moyal equation dAl 8b . As we already saw in paper I, this is not true for the distribution function II e of the 
excursion set formalism, unless one use a sharp filter in momentum space and the theory is gaussian. Already for gaussian theory 
and a different filter, we saw in eq. (83) of paper I that the equation satisfied by II e , besides the Fokker-Planck operator, contains 
complicated non-local terms, coming from boundary terms at the upper integration limit S c . The same happens, of course, when 
we include the non-Gaussianities. 



B. TERM-BY-TERM COMPUTATION OF II (3 < L) 



In SectionlXTIwe showed that 



d5„Tl ( ^\0;6„;S„)=±- jffi 



1 _ £l ) e -^/(2S„) _ 



(Bl) 



Our derivation used the fact that we could replace the sum over i,j, k of djdjdk in eq. d47| i by d 3 /d5l. It is instructive to reproduce 
this result by evaluating separately the various terms in the sum, and using the perturbative formalism of paper I. We then split 
YTi jk=\ m *° tne following terms: (a) i = j = k = n. (b) i < n,j = k = n. (c) i = j <n,k = n. (d) i < j <n,k = n. (e) YlTjk=v eacn one 
with its own combinatorial factor. We denote the corresponding contributions to n (3 L) as n (3,La) , n (3 Lb) , etc. and, for simplicity, 
we also use the notation 

rSc 

I {a) = I d5 n nf™(0;S n ;S n ), (B2) 



and so on. As in the computation of the term proportional to A, ; <9,<9, in paper I, we find that the various contributions in this 
computation can be separately divergent in the continuum limit e — > 0, while their sum is finite, as it clear physically, and as we 
already know from our derivation in Section 13.11 Indeed, the great virtue of the derivation performed in Section 13.11 using the 
trick of replacing the sum over didjdk with derivatives with respect to S c as in eq. d53l l. is that it directly gives the sum over all 
combination of indices, thus providing directly the total finite result, and bypassing all problems of divergences that appear if one 
compute separately the terms corresponding to different combination of indices. 

The separate terms can however be computed using the technique developed in paper I, with the finite part prescription. The 
term has already been computed in eq. ( |59l l, and we have seen that it vanishes. The term (b) is obtained setting j = k = n in 
eq. d47l i. and taking into account a combinatorial factor of three corresponding to the three way of choosing which index, among 
(i,j,k), is not equal to n, so 



n< 3 ' Lb) (<5 ;£„;S„): 
Using eqs. (l65l l and i 



(Si 



n-1 
;=1 



d8 x ... dS n ^ aw*™ = -^lJ2 d n [ n ? m (<5o; $c; sonfU; s n ; s„ - Sd] 



j(b) _ _ 

The term (c) gives 



2tt 



On 



Sc(Sc-S n ) 



exp 



S 2 (8 c -5 n f 



2S t 2(S n -Sd 



J 8„=8, 



2nS 3 J 2 



§ 1 



-sliasn) 



so 



n^Vo;5,,;S,0 = -^i> 



d6 l ...d5 n . 1 d 2 W gm , 



i (c) = --j l 22 ft (nr(<5o;^;5,)nf m (<5 c;( 5„;5„-5,))] 5n=5i 



(B3) 



(B4) 



(B5) 



(B6) 



This expression is analogous to the one that has already been computed in eqs. (B13)-(B16) of paper I, and it is purely divergent 
as 1/ y/e, with no finite part, so J^VU^] = 0. 

The term (d) is slightly more complicated, since it requires the a regularization described in appendix B of paper I for extracting 
the finite part. Setting i < j < n,k = n'm eq. d47T > and taking into account a combinatorial factor of six, we get 



We use eqs. 



and 

j(d) _ _ 



n—2 n— 1 

nf> Lc) (<5 ; 5„;S„) = -{5l)d n n ? m(,5o; &c\Sijnv\8 c \ 5 c ;Sj - Sdll^iS^Sn - Sj) . 

i=\ j=i+l 

together with IT| m (5 c ;(5 c ;5) = e/tV^S 3 / 2 ), see eq. (1 12) of paper I, and we get 

■ >" • ■' Sl A 2 (S,-S„f 



7T V^F S„^8- J J Sj S y\ Sj _ Si )VHS„ - Sj?/ 2 



exp 



25,- 2(Sj-Sd 2(S„-Sj) 



(B7) 



(B8) 



where A 2 = ae regularizes the integral over dSj when Sj ~ * Sf, and we want to extract the finite part as A — > 0. Observe also that 
here one must be careful not to interchange the limit S„ — > S~ (which comes from the fact that the integral over d8„ from — oo to 
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5 C of n^ 3,Lc - ) is performed integrating by parts of the derivative d n that appears in eq. ( IB7t ) with the integrals over dS, and dSj. 
The integrals can be carried out using the identity given in eq. (A5) of paper I, and we get 



lid) . 



2 M> ' 



2tt 



'3/2 V A 



+ 1 exp 



(<5 C +A) 2 
2S n 



(B9) 



This has a part divergent as 1 /A, i.e. as 1 / -y/e, that will combine with the similar divergences from the other terms, a part finite 
as A — * 0, plus terms 0(A) that vanish in the continuum limit. Extracting the finite part we get 



TV\l (d) \=-2 



(ft) 



S„ 



(BIO) 



Finally, the term (e) can be computed with the by now usual trick of replacing YTijk=\ w i tn d 3 /dS 3 , and we get 



ttf Le \6 ;6 n ;S n )=- 



3\ "- 1 r^c 



E 



(d 3 n )4V2 1 



d5i...dS n -idjdjd k W sm = 
(2S c -6„) 2 



(ft) 9 3 



c5/2 



exp 



6 « 

(2<5 c -£„) 2 



and from this we find 



/(e) _ 



(Si) 



3 V2^sl 



3/2 



S 2 C /(2S„) 



(BID 



(B12) 



Summing up eqs. ( IB4b . ( IB 101 ) and ( IB 12b we get back eq. ( IB lb . as it should. 
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